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Towards a Bose glass transition in an optical Penrose quasicrystal 
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Abstract. - We study numerically a 2D Bose-Einstein condensate in a quasiperiodic array of po- 
tential peaks, assumed to be generated by superimposing five blue detuned laser beams. By using 
a Bogoliubov ansatz for the excitations we show that the system approaches a gap less, insulating 
phase upon increasing the potential, consistent with a Bose glass phase. The characteristics of the 
transition in terms of phase correlations, oscillatory modes, and superfluid fraction are discussed. 



Introduction. — The existence of a quasiperiodic 
1 structure was first observed in a rapidly quenched metal 
, alloy by Shechtman et al. [1 . A sharp peak in the Bragg 
scattering proved the existence of long range order. The 
system, however, lacked translational symmetry, under- 
mining the established belief that an ordered structure 
must be periodic. This novel configuration was generalized 
by Levine and Steinhardt [2 , who formalized the concept 
, of quasiperiodic crystals (QC): these formations exhibit 
long range order without having translational symmetry, 
acting both as crystals, admitting a metal-insulator tran- 
sition, and as amorphous solids, with peculiar localization 
properties and fractal spectrum [3 . 

In cold atom research, the investigation of the effects of 
repulsion between bosons when a disordered potential is 
. present is an ongoing quest. In the late 80s it was spec- 
ulated [4j[5] that a new phase of matter exists under this 
condition: The Bose glass (BG). A series of experiments 
in ID [6j[7] showed that such a phase exists also for a 
quasiperiodic potential. This experimental research has 
been accompanied by a corpus of analytical and numer- 
ical articles [8HT4] that helped explaining the properties 
of the glassy phase. In particular, it is predicted that the 
BG should appear in more than one dimension, and ex- 
periments are being done to verify this. 

In this Letter we want to investigate the BG phase in a 
2D QC lattice. The work of Guidoni et al [15j[T6] opened 
up for the possibility of creating an optical QC using an 
appropriate configuration of lasers. In particular, a con- 
figuration with five beams creates a ten-fold symmetric 
structure, similar to the Penrose tiling |17|ll8j. For a non- 
interacting system, this configuration is known to have 




Fig. 1: Penrose potential generated by 5 opposing laser beams, 
as in Eq. (pQ) with Vo = 1 Er. In our configuration, the polar- 
ization of the laser beams is directed towards the reader. 



eigenstates which are neither extended (like Bloch func- 
tions), nor exponentially localized. Surprisingly, few works 
have been done to study the interacting system. Among 
these, a numerical study by Sanchez-Palencia and Santos 
[T9] has shown a that a Penrose-like potential is able to 
inhibit the diffusion of a BEC. We wish to fill the gap in 
the literature and show that the quasiperiodic lattice is 
compatible with a gapless insulating phase, as in the Bose 
glass. 

The optical realization of a QC gives a choice of the 
sign of the potential: the lattice can be a quasiperiodic 
array of wells or peaks. We choose the latter, which can 
be realized with a blue detuned lattice, because in such a 
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Fig. 2: Condensate density (top row) and spatial average of the correlation function as in Eq (0, for g — 10 with Vb = 1, 2, 3, 4 
and 4.5 ^r. In this colorscale blue is the lowest value, red is the highest. 



potential there are no potential wells and a Mott insulating 
phase can not exist. The possible loss of coherence of the 
bosons will then depend solely on localization effects due 
to the quasiperiodicity. 

Penrose tiling. — A quasiperiodic optical lattice can 
be generated in ID by means of superimposing two laser 
beams with incommensurate wavelengths [6]. In two di- 
mensions it is possible to create a tenfold symmetric qua- 
sicrystal structure using five laser beams of equal intensity, 
lying in the same plane with an angle 2 7r/5 between them. 
The resulting potential is given by the expression [20 



V QP (x,y) = V 



E* 



i e*e 



i(ki-r+0i) 



(1) 



where r = (x, y), £{ is the relative dimensionless intensity 
of each laser beam, ki, <pi, and 6{ are the corresponding 
wavevector, phase, and polarization. Vq is an overall in- 
tensity. The wavelength d = 2ir/\k\ of the lattice sets 
the characteristic length of the system. In the following, 
we measure all the energies in terms of the recoil energy 
Er = h/md 2 , and the unit length is d. Moreover, we 
choose £i = 1, (pi = 0, and the various polarization vec- 
tors point all in the same direction, perpendicular to the 
xy plane. In Fig. Q] we plot the resulting potential for 
Vb = 1 Er. This potential appears as a quasiperiodic 
structure of peaks in the plane. 

We consider a quasiperiodic potential combined with a 
harmonic confinement, so that the external potential in 
the xy plane is 



V(x,y) = V QC (x,y) + -uj 2 (x 2 + y 2 ) . 



(2) 



Bogoliubov approach. — The 2D gas of bosons con- 
sidered here will in the absence of a quasiperiodic poten- 
tial form a Bose-Einstein condensate; in the anticipated 



Bose-glass phase, long-range coherence will be lost so that 
the gas can be considered a quasicondensate [21] . In both 
cases, the system is accurately described by a (quasi-) con- 
densate wavefunction ipo(x,y), accompanied by quadratic 
fluctuations described by Bogoliubov modes [22-24 . The 
(quasi-) condensate wavefunction is governed by the well- 
known 2D Gross-Pitaevskii equation (GP2D) 



:V 2 ^o + 9 |^o| 2 ^o + V i/j = /i ^o . 



(3) 



In this 2D approximation, the wavefunction is supposed 
to be constant in the direction perpendicular to the xy 
plane, extending for a length equal to a±. The value of 
g = g( 3D )/ax is the scaled interaction strength among 
the particles in 2D, V(x, y) is the external potential, and 
\i is the chemical potential. The density of the (quasi-) 
condensate is given by n = |^o| 2 - Up to second order in 
the many-body Hamiltonian the atoms outside the (quasi- 
) condensate occupy excited states which are determined 
by solving the Bogoliubov equations 



1 2 

- - V + Vj + 2 n g — fi 1 Uj — ngvj 

- - V 2 + Vj + 2 n g — fi ) Vj — nguj — —ujj Vj .(4) 



where, as usual, the normalization J Ui Uj — V{ Vj = Sij 
is enforced. As was found in Ref. [14j[22], the one-body 
correlation function can be expressed in terms of the Bo- 
goliubov excitations as 
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Fig. 3: Logarithm of the rotation average of the correlation 
function, as in Eq. (|5J). 



provided ipo is real. In this formula, ijj is the 
many-body boson operator, and for brevity we write 
Uj,Vj,n and u'^v'^n' instead of Uj(r),Vj(r),n(r) and 
Uj(r),v'j(r), n'{v). Nj is the occupation number of the j th 
excited state, as determined by the Bose distribution. In 
current experiments, this occupation can be made so small 
that it can be neglected with respect to the depletion from 
the (quasi-) condensate given by the interaction. In the 
following we work with Nj = (as for vanishing temper- 
ature) and the only contribution to the correlation is the 
first term in the sum. 

We remark that the Bogoliubov approach is a gapless 
ansatz for the excitations in a Bose system. This approach 
works for the gapless Bose glass phase, but it fails to detect 
the Mott insulator transition because the latter phase has 
a gap. The choice of our potential inhibits the formation of 
a Mott phase, and the only insulating phase in the system 
is provided by the quasiperiodic pattern of the lattice. 

We find the ground state of the GP2D by using an 
imaginary time evolution with the Fourier split opera- 
tor method. Subsequently, we diagonalize the Bogoli- 
ubov equations by employing the ARPACK libraries with 
Tchebychev polynomial acceleration [25] . We choose a sys- 
tem size of L = 80 <i, and the square grid has 480 points 
per side. The matrix to diagonalize has therefore the di- 
mensions of 2 x 480 2 ~ 4.6 x 10 5 points per side. This 
matrix is not banded, since we have chosen to represent 
the Laplacian using the Fourier transform. Therefore, in 
order to have convergence in reasonable time, even with 
the polynomial acceleration, we must make a cutoff at 100 
eigenvectors. We shall see that this cutoff is adequate. 

We consider a system with g = 10; this relatively weak 
interaction can be realized with, e.g., a gas of ~ 10 5 atoms, 
confined in a transverse dimension a± ~ lOOnm, with a 
scattering length of ~ 10 _4 nm, which we propose could 
be realized using Feshbach resonances. The harmonic po- 
tential strength in the xy plane is u = 8 x 10 -2 Er. The 



reason for these parameter values is numerical limitations; 
in principle, there is nothing to prevent the BG transition 
to occur also for larger values of g. Also, the numerical 
results are limited to Vo < 5.75 Er: Beyond this param- 
eter regime, an unattainably high spatial resolution was 
required in order to correctly represent the mode func- 
tions. 

The top row of Fig. Oplots the density of the condensate 
while changing the strength of the lattice from Vo = 1 £?r 
to 4.5 ^r. In the bottom row we plot instead the spatial 
average of the correlation function [26] 



5 (r) = i|dR<V>t(R + r )</;(R)). 



(6) 



This quantity is directly related to what is imaged in ex- 
periments; the density after time of flight is to a good 
approximation equal to the momentum space density pk, 



Pk 



V 



drdKe lk - r g(r). 



(7) 



As we can see, for low values of Vo the characteristic width 
of g(r) is the entire condensate, indistinguishable from the 
average correlation function of a pure condensate in a har- 
monic potential. When Vo ~ 3£^ this correlation func- 
tion starts to shrink, and eventually decays rapidly when 
getting farther from the center. 

This behavior is particularly evident when we plot the 
angular average of g(r) 



g(r) 



/ 



dB_ 

2tt 



sM), 



(8) 



as is shown in Fig. [3] When a weak lattice is present, 
the logarithm of the correlation function decays on a scale 
comparable with the total size of the system. Upon in- 
creasing the lattice, a central peak appears where the log- 
arithm of the correlation decays linearly, on a scale shorter 
than the system size. This is one of the signatures of a 
crossover to an insulating phase, as reported by Deissler 
et al 0. 

Elementary excitations. — It is known that the 
Bose glass phase is caused by low-lying excitations that 
can flip the phase of the quasicondensate with a little 
amount of energy. One of the signatures of the glassy 
phase is the progressive lowering of the excitation spec- 
trum. 

In Fig. H] we plot the energies of the lowest dipole and 
quadrupole excitations. We see that the raising of the 
quasiperiodic lattice separates each of the modes into two 
branches. More importantly, the energy of the excitations 
drops significantly. These modes can be excited by im- 
posing dipolar and quadrupolar deformations in the con- 
fining potential. In particular, a strain can be given to 
the quasicondensate by a slight change in the harmonic 
confinement, 



-u (x + y ) 



i(c, + e) 2 * 2 + i(c,-e)V. (9) 
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Fig. 4: Dipolar and quadrupolar excitations. Right panel: The 
four lowest excitation energies, in units of Er = h 2 /md 2 , as 
functions of the lattice strength Vq. Lower left: Bogoliubov 
amplitude v± for the first excited mode; upper left: Bogoliubov 
amplitude ^3 for the third excited mode, both computed for the 
case Vq = 1 ^r. 



Such a deformation will predominantly set off an exci- 
tation with an irrotational velocity field. For Vo ~ 0, 
the characteristic frequency of this mode is seen to be 
Co « 0.13 £?# « 1.6cj. This is the scissor mode, described 
in Ref. [27]. Using a hydrodynamical approach, Ref. [27] 
shows that a superfluid admits quadrupolar excitations 
(scissor modes) in a harmonic potential, and the lowest 
frequency for this type of excitation is in the strongly in- 
teracting limit equal to \^2uj; for our finite value of g, it is 
slightly higher at Vo = 0. Another important result of Ref. 
[27] is that a normal fluid would dampen this excitation 
in a finite time. We argue that the disappearance of the 
scissor mode provides a way to determine the phase transi- 
tion to the Bose glass, since when the system is completely 
normal the quadrupolar excitations are short-lived. 

In order to compute the lifetime of the scissor modes one 
should compute the imaginary part of the stress-tensor- 
stress-tensor response function. This cannot be done in 
the simple Bogoliubov approximation, because even for 
a uniform system the response function involves an inte- 
gral that is ultraviolet divergent. In order to lay down a 
microscopic theory of the quadrupolar excitations a renor- 
malization approach seems compulsory, but such a work 
is outside the scope of the current letter. 

Superfluid fraction. — The clearest signature of the 
Bose glass is the fact that the superfluid fraction of the 
density vanishes, and the system acts as a normal fluid. 
The superfluid fraction is the subset of the fluid that can 
move only according an irrotational velocity field, while 
the rest of the gas constitutes the normal fraction. In- 
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Fig. 5: Normal fraction of the system as function of Vo (in 
units of Er = tf/md 2 ), according to Eq. (|14|) . The circles 
correspond to the ratio I /hot computed according to Eq. (|14[) . 
In order to estimate where the superfluid part disappears, we 
plot a spline interpolation of the whole data, a linear fit of the 
last two points, and an exponential fit of the last three points. 



spired by the work in Ref. [28], we compute the normal 
fraction by looking at the response to a rotation of the 
system. The external potential can be rotated at a small 
angular velocity Q. In the rotating frame the Hamiltonian 
of the system appears as 



H — H — \ I L z , 



(10) 



where the angular momentum operator is (using hats on 
second-quantized operators but not on first-quantized op- 
erators) 



L z (r,t) = L z (r) ft(r',t')$(r,t) 



r' —>r ,t' -tt 



(ii) 



and r = (#, y), so that L z (r) = — ih(xd y — yd x ). 

If the stirring is slow enough, superfluid vortices can- 
not be generated, and only the normal part takes part in 
the rotation of the fluid. In the limit of vanishing Q, the 
increase in energy due to the rotation is given by 



Eq — Eq 



1 



Intf 



(12) 



where I n is the moment of inertia of the normal part. 
Using standard perturbation theory, the normal moment 
of inertia is given by second order perturbation theory. 
Expressed as integrals over the imaginary time r = —it, 
it reads 



I n = -5 I dn dr 2 



dr 1 dr 2 (L z (r 1 ,T 1 )L z (r 2 ,T2)) , (13) 



with /3 = 1/T, the inverse of the temperature. Notice that 
this expression is the q = 0,cj = limit of the angular- 
momentum-angular-momentum response function, as de- 
scribed by Pines and Nozieres [28 . By applying Wick's 
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theorem we obtain, in the Bogoliubov approximation 

^ = 2j2^^ + Y,-> ( 14 ) 



' UJi + Uj ^ UJi 



with 



*#j 



Xij = ( ViL z Uj) ( / Uj L z Vi) 



+ ( / ViL z Uj)( / uiL z Vj) 



(15) 



where the last integrals are only over the spatial coordi- 
nates, since the (imaginary) time has been integrated out 
in the evaluation of the Matsubara frequencies [14 . The 
normal fraction is given by the ratio I n /I tot , where I t ot 
is the total momentum of inertia. This ratio is plotted in 
Fig. [5] Note that the sum over the excitations is limited 
by the cutoff in the diagonalization procedure. On the 
other hand, we checked that this restriction of the num- 
ber of modes does not alter the result: The highest 10% of 
the states included were seen to contribute less than 1% 
of the sum. This is expected, since the Bose glass phase 
depends on the lowest lying excited modes. 

According to our analysis, the normal part starts to 
increase when the correlation function shows an exponen- 
tial decay on a scale of the whole system (Vb ~ 5^r). 
We cannot see the full transition because no numerical re- 
sults could not be obtained for Vb > 5.75 £?r, as discussed 
above. However, an assortment of extrapolation methods 
- Fig. [5] shows results of spline, exponential, and linear 
extrapolation - indicate that the superfluid part should 
vanish between Vo = 6 Er and Vb = 9 Er. Moreover, note 
that the transition to the glassy phase is not sharp, but 
it appears as a smooth crossover from superfluid to insu- 
lator. This effect is due to the finite dimensions of the 
system, and for the ID Bose glass it has been experimen- 
tally observed in Ref. [7 . 

Conclusions. — We have seen evidence suggesting 
that a a blue detuned optical quasicrystal generates a 
phase transition in a Bose gas, from a superfluid to an 
insulating phase. The signatures of this transition are an 
exponential decay of the correlation function and an in- 
crease in the normal part of the gas. Since we have chosen 
a gapless ansatz for the excitations, and a the blue detuned 
lattice without potential wells, we know that this phase is 
not a Mott insulator. The transition is due to the specific 
shape of the potential and the localization effects it causes 
on the Bose gas. Summing up, the quasiperiodic pattern 
of the lattice leads to a normal, gapless state, compatible 
with the description of the Bose glass. 
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